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Abstract 

We propose a novel method for the extraction of unresolved point sources from CMB 
maps. This method is based on the analysis of the phase distribution of the Fourier 
components for the observed signal and unlike most other methods of denoising does not 
require any significant assumptions about the expected CMB signal. The aim of our paper 
is to show how, using our algorithm, the contribution from point sources can be separated 
from the resulting signal on all scales. We believe that this technique is potentially a very 
powerful tool for extracting this type of noise from future high resolution maps. 

Subject headings: cosmic microwave background, cosmology, statistics, observations. 

1 Introduction 

Observations of the Cosmic Microwave Background (CMB) is fundamental for our under- 
standing the primordial inhomogeneity of the Universe. After the successful COBE experi- 
ment, attention has been focused on the investigation of small scale perturbations, that can 
provide unique information about the most important cosmological parameters. One of the 
major problems in the modern CMB cosmology is to separate noise of various origins (such 
as dust emission, synchrotron radiation and unresolved point sources (see e.g. Banday et 
al. 1996)) from the original cosmological signal. Many authors have already applied various 
methods such as Wiener hltering (Tegmark and Efstathiou 1996, Bouchet and Gispert 1999), 
maximum entropy technique (Hobson et al. 1999), radical compression (Bond et al. 1998), 
power filtering (Gorski et. al. 1997, Naselsky et. al. 1999) and wavelet techniques ( e.g. 
Sanz et al. 1999) to extract noise from the CMB data. 

All of these techniques have been tested for removing the noise from the real observational 
data. It is necessary to note that, for different strategies and for different experiments, 
different schemes could be chosen as most appropriate. The choice of the algorithm also 
depends on the particular type of foreground emission to be extracted. 
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The aim of our paper is to overcome the problem of detecting and extracting the back- 
ground of unresolved point sources from the original map. The measured signal in the real 
observational data is always smoothed with some filtering angle 9f because of the final antenna 
beam resolution. Therefore, unresolved point sources could make a significant contribution 
to the resulting signal on all scales. This type of noise should be removed from the original 
map before any subsequent analysis is made. 

Recently (Cayon et al. 1999) have proposed the use of isotropic wavelets for removing 
noise in the form of point sources. Their technique is based on the fact, that the field in 
the vicinity of the source should be in the form of the antenna profile. Unfortunately the 
Gaussian CMB field can also form real peaks with the same profile, so that a lot of 'artificial 
sources' could be found using this technique. Besides, the antenna profile is not necessarily 
isotropic (indeed, as a rule it is very anisotropic). Therefore, isotropic wavelets should not 
be considered as an absolute cure against such a type of noise. 

In this paper we consider an approach, which is based on the distribution of phases. 
The idea of using phases of random fields was introduced by A.Melott et al (1991); Coles 
and Chiang (2000a, b) for the Large Scale Structure formation in the Universe. Below we 
develop the phase-amplitude analysis method for investigation of the CMB anisotropy and 
foreground. The outline of the paper is as follows. In section 2 we briefly review the basic 
definitions, consider a simulated one-dimensional scan of the CMB first with a single point 
source, then with a background of such a sources. In section 3 we generalize our results 
into two-dimensional maps. Finally, we suggest an algorithm for denoising. In section 4 we 
discuss the results and potential of the method for analyzing high resolution maps. 

2 Point sources in one-dimensional scans. 

In this section we consider ID CMB scans with a background of point sources. This approach 
could be very Tiscful for data analysis of one-dimensional experiments with high resolution 
(such as RATAN 600). We extend this discussion to two-dimensional experiments (such as 
the new generation of interferometer experiments) in section 3. The investigation of point 
sources is especially easy in one dimension, can be easily generalized into two-dimensional 
maps and will help us to understand the advantage of the proposed technique. 

Definitions 

In ID the deviation of the temperature from its mean value AT = T — (T) in a scan is 
described by the simple Fourier series: 

Ar(^) = 5Z cos(A;^) + bk sm{ke) (1) 

k 

where k is an integer number and 9 can be expressed in terms of of the real angle on the sky 
{Osky) as follows: 9 = ^^6sky Here Otot means the total length of the scan. 

The detected temperature fluctuations AT can as usual be naturally divided into two 
parts: cosmological signal and noise: 
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Ar(^) = AT' {9) + Ar"(0) (2) 

where s and n denote signal and noise respectively. Therefore, the Fourier transform compo- 
nents Ofe, bk can be also expressed as a sum of Fourier decomposition of these two terms: 

h = bi + hl. 

The statistically isotropic distribution of the CMB temperature anisotropy is supposed 
to be in the form of a random Gaussian field with the power spectrum Pc'MB{k), which 
determined by the appropriate cosmological model. The coefficients a|,5| depend on the 
spectrum of the CMB, the antenna filtering function F{6 — 0*,0f) and the actual real- 
ization of the random Gaussian process on the sky. In general, they obey the formulae: 
(a|a|,/) = = 6kk'F{k, kf)PcMBik)- Here, F{k,kf) is the Fourier transform of the fil- 

tering function and ^/ is the antenna resolution angle, fc/ is a wavenumber which corresponds 
to this resolution: kf = l/9f. In our simulations we use the usual expression for ak,bk: 

1 1 
af^ = akF2{k,kf)P^^j^{k), 

bl = PkFHk,kf)PSMBik), 

where ak,(3k are independent Gaussian numbers with zero mean and unit dispersion. 

In this paper we consider the noise in the form of isolated unresolved point sources. 
This means that the average distance between sources is larger than the resolution scale 9f. 
Therefore, the shape of the 'noise' field around the point source determined by the filtering 
function F: 



ATn{9) = / J2^jd{9*-9)F-2i9*-9,9f)d9* 



= Y^.^^F-2{9-9,,9f) (5) 

where 7j, 9j are the amplitude and the position of the j-th point source, respectively, and Npg 
is the total number of point sources in the considered scan. According to equation [5], the 
Fourier components of the noise can be described by the following very simple and convenient 
formulae: 



fifc = E -fjCos{k9j)F2{k,kf), 
bl= j:ijSmik9,)F-2{k,kf). 

For further investigation we have to introduce the phase: ipk of the k-th harmonic. Using 
equations [1,3] one can write: 
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= arctan 



ak 



arctan 



H + K 



(7) 



If the resulting field at the scales k is dominated by the Gaussian CMB signal {Sk/Nj- >> 
1), then ifi^ ~ arctg{b1^/ a1?). In this case the phases of the k-th harmonics are random 
independent uncorrelated values, uniformly distributed from to 27r. On the other hand, 
if the signal at these scales is much smaller than the noise, then the distribution of phases 
is determined by the positions and amplitudes of point sources on the scan. In Fig.l, we 
present the spectrum of CMB in one dimension PcMsik) for the standard CDM model 
together with the spectrum of point sources. Both spectra are smoothed with the Gaussian 
filtering function F(k,kf) = exp(— ^^). It is well known that the CMB signal disappears 

when k becomes larger than some value kd- This value corresponds to the damping scale 
of the CMB fluctuations. Therefore, at the small scales the resulting field is dominated by 
the noise. Note that kd should not be necessary interpreted as the damping scale. Roughly 
speaking, this is the scale where noise from sources becomes larger then the CMB signal. 




Fig. 1 The power spectrum of the CMB (solid line) for one-dimensional 
(360°) scan together with the spectrum of point sources (dashed line). 



It is easy to see, from equations [4,6,7], that the process of smoothing does not change 
the phases of the primordial signal. The filtering function F{k,kf) has simply disappeared 
from the right hand side of the equation [7]. Therefore, if kf > kd, we have the possibility 
of measuring the phases only for high k values of the noise. Below we describe how the 
information about the phase distribution for high values of k can be used for very precise 
detection and extraction of the contribution from the sources for all values of k 

Detection of a single point source 

Let us consider the simplest example by dealing with a single unresolved point source on 
the scan. In order to remove the contribution from this source, we have to know its precise 
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location 9i and amplitude 71 (see Fig.2). 
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Fig. 2 Upper panel: simulated CMB field on 100° scan (1° corresponds 
to w 0.03° on the sky) (dashed line) and the field from a single point source 
(solid line). Lower panel: the same as the upper one, but with better 
resolution. The field in the vicinity of the point source behaves like an 
ordinary Gaussian fluctuation. 



The contribution from this source to the resulting field according to equation [5] is then: 

k=l 

where kmax is the maximum value of k that can be detected in the experiment. 

As has been already mentioned, for k larger then some value kd, phases ipk are just the 
phases of the point source. From equations [6,7], we obtain: 

ifk = mod2Tr{k9i) (9) 
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It suffices to have only two phases (for example ^pk and ^k+i-, k > kd) to find the location 
of the source ^i: 

Oi = ^k+i - (10) 

In Fig. 3 we show the behavior of the phases <pk, 1 < k < kmax together with the phases of 
the source. For small values of A:: k « k^ the phases are distributed uniformly and at large 
k we can definitely see the regular structure that is consistent with equation [9]. 
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Fig. 3 Tlie phases of the point source (circles) and phases of the resulting 
signal: CMB + Point source (crosses). 



In Fig. 4 we also show the positions of maxima for all harmonics. Location of the maxima 
for the k-th harmonic can be found by the formulae: 

nk ^ (^fc + 27r * n 

where n is an integer number. The straight vertical line points to the location of the source 
because one of the maxima in each harmonic is coincident with this location. 
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Fig. 4 Positions of maxima for each harmonic. Each point represents the 
positions of maxima for the k-th harmonic. For smaU k they are distributed 
uniformly (according to the Gaussian distribution of the CMB). The large 
dot shows the location of the source. 



The remaining part of the problem is to find the amplitude - 71. Let us defined the field 
AT^'^{0) as a part of the field AT{9) that consists only of the high harmonics: 

AT^d{e) =^ak cos{k0) + bk sm{k9) (12) 

k=kd 

Using the formulae [8], we now can write down the obvious relation: 

^^ = AT^^{ei)/ J2 Fik,kf) (13) 

Therefore, according to [3,6], we have found the contribution from this source to all 
harmonics from A; = 1 to A; = kmax- 

Background of point sources 

In this subsection we generalize our algorithm to the case where there are an unknown 
number of point sources in the considered scan. In a situation like this, we have to find not 
only positions and amplitudes of each source but also the total number of them: Npg. 
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We believe that many different techniques based on the results of the previous subsection 
could be proposed to solve this problem. We suggest a simple iteration scheme. As has been 
already noticed above, we can consider the field AT^''-, which consists only of high harmonics. 
Therefore, only point sources make a contribution to this field: 

AT''d{e)= J2 a]^cos{ke) + b]^sm{ke) = 

N , (14) 

= E7j E F2{k,kf)cos{k{e-ej)) 

j=l k=kd 

We now introduce the filter function L{k) = F{k — kd, ki)/F{k, kf) and consider the filtered 
field: 

AT^d ^ ATj^'^L^k) (15) 

According to [14,15], one can write: 

Af''d{e) = 

Nps kmax~k^ ^ 

Eijcosikdie-ej)) E cos{k{e -ej))F{k,ki)2 , . 

j=i k=i y^'^) 

kmax—k^ , 

- sm{kd{e - %)) E sin(fe(e - ej))F{k, ki)2 
k=l 

If we can put ki « kd « kmax, then the second term on the right hand side of equation 
[16] is small and: 

Nps kmax—k^ 

Af^^{e)^Y.^^ E cos(fc(^-0,))F(fc,fcO^ (17) 

3=1 k=l 

This equation is very close to [8] and, therefore, the procedure of filtering gives us the possi- 
bility of 'localizing' the field in the vicinity of a point source. 

In reality equation [17] is not quite correct because kd and kmax are values of approxi- 
mately the same order and, therefore, peaks, that are more or less close to each other can 
interfere (fig. 4). This is the reason why we choose the iteration technique to remove point 
sources. 

We propose the following algorithm. Let us construct the field ATJ = AT^'^ and find its 
highest maximum. This maximum most probably corresponds to the most powerful isolated 
point source on the scan. The position and value of this maximum give us the location 9\ and 
the amplitude 71 (eq[13]) of this source. After that, we construct the field ATf'' 'without' 
this point source: 

kmax 

AT^(0) = Af^^e) - 71 E cos(A:(e - ^i))L(A;, fc;)^ • (18) 

k=kd 

The contribution from this source to the field and its interference with other sources is now 
removed. This allows us to find more precisely the next highest maximum. Therefore, we 
apply the same procedure to the field Arf"* and find ^2, 72 and so on (Fig.5). 
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Fig. 5 The iteration sclieme. Eacli panel represents the residuals AT^'^ 
from the initial field AT^-^ = AT*'' after the i-th iteration. 



We perform these iterations until the dispersion a'l'^ {{<^i^)^ = ii^T^'')'^}) becomes signif- 
icantly smaller then a^'' (Fig. 6). The total number of iterations that is needed to significantly 
reduce the initial dispersion gives us approximately the number of point sources Npg and each 
iteration gives the location 9i and the amplitude 7^ of the i-th source. Note, that 



(19) 



and roughly speaking, in Fig. 7 we can see the cumulative distribution of point sources over 
the power 7. Finally, since we have the position 6i and amplitude 7^, the contribution to the 
field from all point sources may be removed in the same manner, as was done for a single 
point source in the previous subsection. 
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Fig. 6 The decrease of the dispersion for AT'''^ with each iteration. 
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Fig. 7 The final result. The initial field of point sources (solid line), 
restored field by our method (dashed fine), and residuals (dotted line). 
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3 Point sources in two dimensions. 



In this section we briefly describe our results in two dimensions. Without loss of generality 
we may consider a small region of the sky and assume that the geometry is approximately 
flat. Under this assumption, the part of the detected signal which is determined by the noise 
associated with Npg point sources can be represented according to the previous section by 
writing: 

ATn{x) = E Ofc cos{kx) + bl sin(fe) = 

N.S ^ (20) 

E ijY.F-^{k,\kf\)cos{k{x-Xj)) 
J=i k 

where Xj is the position of the j-th point source in the Cartesian coordinate system and \kf\ 
corresponds to the antenna resolution. Analogously to the one-dimcnsional case, this field 
should be filtered with some appropriate function. The convenient filter function L(k) that 
we use in this case is as follows: 



L{k) = F-^{{k), \kf\) if \kd\ < \k\ < \kmax\, /21) 

L{k)=0 if \k\ < \kd\ or \k\ > \kmax\. 

According to [15,20,21] one can write: 

^ps I ^ I ~ I ^max I 

Af'='^(^) = 5:7, E cos(fc(f-x,-)) (22) 

l^l=|fedl 

Therefore, the filtered function /S.T^'^{k) at the point x = Xj has a peak with amplitude equal 
to the power of j-th point source times the number of modes that we can use for data analysis 
in the appropriate experiment. 

In our simulations of the signal+noise we use the standard CDM model and background 
of 100 point sources randomly distributed over the 10° x 10" map. Without loss of generality 
we use the simple symmetric Gaussian antenna profile. All these calculations, of course, could 
be done for any arbitrary antenna beam. In Fig. 8 we show the map of the CMB together with 
the maps of noise, CMB+noise and the filtered map of CMB+noise. The last one shows us 
more or less clearly the positions and powers of point sources. The significant anisotropy that 
appears in the last map occurs for the following reason. According to formulae [21] we use 
only the set of harmonics fei, /c2, that obey the relation fe^ + fe^ > fc^. Therefore, the number 
of horizontal and vertical waves is larger than the number of waves in any other direction. 
(This problem does not occur if we use spherical harmonics with I > Id). 

We apply the same iteration technique as in the previous section and therefore separate 
noise from the cosmological signal (Fig.9). It is necessary to note that each iteration removes 
an appropriate point source at the beginning of this process for the most powerful and sep- 
arated sources. For the weaker sources, additional iterations are needed. The signal from 
the i-th source decreases as ~ --^ — ^ , where dr is the distance from the peak (in one 

dimension this dependence is linear). 



12 




CMB + noise CMB + noise 



after filtration 

Fig. 8 Simulated sky maps of 10° x 10°. 
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Fig. 9 Noise maps (i.e. removed sources) after different numbers of iterations. The size and shading of each source is 

proportional to its amphtude. 
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This affects tlie neigliboring peaks and can change their amphtudes. Therefore, this 
approximation works if {{kmax — kdjSvij)'^ « where rij is the separation between the i-th 
and j-th peaks (the i-th peak is the closest to the j-th peak). Otherwise, the amphtudes that 
have been found in each iteration would not correspond to the powers of the sources and we 
therefore have to perform a number of iterations that is larger than the number of sources. 

4 Conclusions 

In this paper we present a powerful method for extraction of unresolved point sources from 
future high resolution CMB maps (such as MAP, Planck, VSA, CBI, DASY, AMI and RATAN 
600). Our method is based on the distribution of phases. The most important advantage 
of our technique is that we do not make any strong assumptions about the expected CMB 
signal as well as about the antenna profile. Most other techniques use the estimated power 
spectrum of the CMB and noise before the data analysis is implemented (e.g. Wiener filter) 
or they require special assumptions about the antenna profile (e.g., wavelets techniques). It 
is worth stressing that, for example, assumptions about the CMB power spectrum can lead 
to incorrect interpretations of the observational data. Roughly speaking, by making such 
assumptions, one runs the risk of generating the result one wants and any discrepancies are 
consider to be errors. 

Our algorithm is numerically very efficient. It is a linear algorithm and requires Nln(N) x 
Nps operations, where N is the number of pixels. Therefore it can be easily applied to the 
analysis of large data sets. We have demonstrated the accuracy which can be achieved using 
our algorithm to remove the contribution from point sources on all scales. We believe that 
this technique is potentially a very powerful tool for extracting this type of noise from future 
high resolution maps. 
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